#Note: this code estimates Models 5-8 as reported in Table 1 of the Research Note.
#   To do so, it reads-in the Stata dta files "MainMalariaDatasetStata11.dta" and "MonMalariaDatasetStata11.dta"
#   and then prints the relevant model output to screen.


################################################################
#############################Set-up#############################
################################################################

#load necessary packages
library(plm)
library(doBy)
library(plyr)
library(readstata13)

#set working directory
setwd("")

#################################
###Read in Data for Models 5-6###
#################################

#import dyadic dataset
dyad.dat <- read.dta13("MainMalariaDatasetStata11.dta")

#############
###Model 5###
#############

##Keep only variables of interest
base.dyad <- dyad.dat[,c("dr_1_at_2_m","Lmalaria_interp_2_m","ddyad","year5")]   

#Convert data to pgmm formal
pgmm.base.dyad <- pdata.frame(base.dyad, index=c("ddyad", "year5"))

#set seed
set.seed(50)

#run and summaize model
gmm.base.dyad <- pgmm(dr_1_at_2_m ~ Lmalaria_interp_2_m| 
                        lag(dr_1_at_2_m , 2:5)+lag(Lmalaria_interp_2_m,2:5),
                      data = pgmm.base.dyad, 
                      effect = "twoway", transformation= "ld", model = "onestep")
summary(gmm.base.dyad)

##Number of observations
base.dyad.na <- na.omit(dyad.dat[,c("dr_1_at_2_m","Lmalaria_interp_2_m","ddyad","year5")])   
length(base.dyad.na$dr_1_at_2_m)

#############
###Model 6###
#############

##Keep only variables of interest
full.dyad <- dyad.dat[,c("dr_1_at_2_m","Lmalaria_interp_2_m","Llog_trade_m","Ldr_2_at_1_m","Llog_gdppc_1_m","Llog_gdppc_2_m","Ls2_m","Lln_cinc_1_m","Lln_cinc_2_m","Ldem_int_m","Ldemocracy_1_m","Ldemocracy_2_m","ddyad","year5")]

#Convert data to pgmm formal
pgmm.full.dyad <- pdata.frame(full.dyad, index=c("ddyad", "year5"))

#set seed
set.seed(50)

#run and summaize model
gmm.full.dyad <- pgmm(dr_1_at_2_m ~ Lmalaria_interp_2_m+
                        Llog_trade_m+Ldr_2_at_1_m+Llog_gdppc_1_m+Llog_gdppc_2_m+
                        Ls2_m+Lln_cinc_1_m+Lln_cinc_2_m+Ldem_int_m+
                        Ldemocracy_1_m+Ldemocracy_2_m| 
                        lag(dr_1_at_2_m , 2:5)+lag(Lmalaria_interp_2_m,2:5),
                      data = pgmm.full.dyad, 
                      effect = "twoway", transformation= "ld", model = "onestep")
summary(gmm.full.dyad)

##Number of observations
full.dyad.na <- na.omit(dyad.dat[,c("dr_1_at_2_m","Lmalaria_interp_2_m","Llog_trade_m","Ldr_2_at_1_m","Llog_gdppc_1_m","Llog_gdppc_2_m","Ls2_m","Lln_cinc_1_m","Lln_cinc_2_m","Ldem_int_m","Ldemocracy_1_m","Ldemocracy_2_m","ddyad","year5")])
length(full.dyad.na$dr_1_at_2_m)


#################################
###Read in Data for Models 7-8###
#################################

#Import monadic dataset
mon.dat <- read.dta13("MonMalariaDatasetStata11.dta")

#############
###Model 7###
#############

##Keep only variables of interest
base.mon <- mon.dat[,c("dr_1_at_2_m","Lmalaria_interp_2_m","cowcode_2","year5")]   

#Convert data to pgmm formal
pgmm.base.mon <- pdata.frame(base.mon, index=c("cowcode_2", "year5"))

#set seed
set.seed(50)

#run and summaize model
gmm.base.mon <- pgmm(dr_1_at_2_m ~ Lmalaria_interp_2_m| 
                       lag(dr_1_at_2_m , 2:5)+lag(Lmalaria_interp_2_m,2:5),
                     data = pgmm.base.mon, 
                     effect = "twoway", transformation= "ld", model = "onestep")
summary(gmm.base.mon)

##Number of observations
base.mon.na <- na.omit(mon.dat[,c("dr_1_at_2_m","Lmalaria_interp_2_m","cowcode_2","year5")])   
length(base.mon.na$dr_1_at_2_m)

#############
###Model 8###
#############

##Keep only variables of interest
full.mon <- mon.dat[,c("dr_1_at_2_m","Lmalaria_interp_2_m","Llog_trade_m","Llog_gdppc_2_m","Lln_cinc_2_m","Ldemocracy_2_m","cowcode_2","year5")]

#Convert data to pgmm formal
pgmm.full.mon <- pdata.frame(full.mon, index=c("cowcode_2", "year5"))

#set seed
set.seed(50)

#run and summaize model
gmm.full.mon <- pgmm(dr_1_at_2_m~Lmalaria_interp_2_m+Llog_trade_m+Llog_gdppc_2_m+Lln_cinc_2_m+
                       Ldemocracy_2_m| 
                       lag(dr_1_at_2_m, 2:5)+lag(Lmalaria_interp_2_m,2:5),
                     data = pgmm.full.mon, 
                     effect = "twoway", transformation= "ld", model = "onestep")
summary(gmm.full.mon)

##Number of observations
full.mon.na <- na.omit(mon.dat[,c("dr_1_at_2_m","Lmalaria_interp_2_m","Llog_trade_m","Llog_gdppc_2_m","Lln_cinc_2_m","Ldemocracy_2_m","cowcode_2","year5")])
length(full.mon.na$dr_1_at_2_m)

